function [ W1, H1 ] = TNMF( duomenys, R, niter)

    M = size(duomenys, 1);
    N = size(duomenys, 2);
    
    Y = duomenys;
    W1 = rand(M, R);
    H1 = rand(R, N);
    
    ee = 1e-6;
    
    for k=1:niter
        
        W1 = ApplyNMFNormalization(W1);
        
        W1 = W1 .* (Y*H1' )./((W1*H1)*H1' + ee);
        H1 = H1 .* (W1'*Y )./(W1'*W1*H1 + ee);
        
        
        W1 = Limiter(W1, R, [0.0 0.65; 0.2 0.8; 0.4 1], ee, 1);
        %H1 = Limiter(H1, R, [0 0; 0.0 0.5; 0 0], 0.5, 0);
        
        %nimble = [nimble dist];
    end
    
    dist = DivergDist(duomenys, W1 * H1);
        
    if dist > 5,
        disp('Not close solution - retrying');
        [W1 H1] = TNMF(duomenys, R, niter);
    end
    
    
end
